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Charmonia are flavour singlet mesons and thus in principle contributions from disconnected quark 
line diagrams might affect their masses, either directly or via mixing with other flavour singlet 
channels. We present a first study that takes both effects into account. We employ improved 
stochastic all-to-all propagator techniques (including new methods) to calculate the diagrams that 
appear within the mixing matrix between the T]' and the 77c. The runs are initially performed 
on Nf = 2 16^ X 32 configurations with the non-perturbatively improved Sheikholeslami- Wilson 
action, both for valence and sea quarks. 
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1. Introduction 

Lattice calculations of the charmonium S'-wave hyperfine splitting tend to underestimate the 
experimental value of 117 MeV [|l|-^. However, in none of these studies all systematics have 
been addressed as yet. One such source of error is the difference of the running of the QCD 
coupling between low momenta ~ mv^ that determine the spin averaged splittings and the high 
momenta > mv responsible for the fine structure. Too few or too heavy sea quark flavours will 
result in a comparably smaller coupling at high momenta, resulting in such an underestimation |^. 
Furthermore, the finestructure is particularly sensitive to short distance physics, necessitating a 
careful extrapolation to the continuum limit. 

Here we will investigate the effect of quark-antiquark annihilation diagrams, the neglection of 
which might be an additional source of the discrepancy. Although heavy, compared to the chiral 
symmetry breaking scale, the rjc state might still sense some remnant axial U{1) effect that could 
result in an upward mass shift, perhaps of a few MeV. Another aspect is the mixing with other 
flavour singlet states, e.g. the pseudoscalar glueball or the rj' meson and its radial excitations. 

We study both effects: quark annihilation diagrams [[l]-^ and mixing with the rj', hence called 
77 since we have Nj = 2. To investigate this mixing we construct a correlation matrix containing 
both light and charm quark operators, the calculation of which requires all-to-all-propagators. The 
eigenvalues and appropriately normalized eigenvectors of this matrix are indicative of the magni- 
tude of mixing between the two sectors. 

2. Simulation details 

Our runs are performed on Nf = 2 16^ x 32 configurations generated by the QCDSF collab- 
oration [||], with a lattice spacing a « 0. 1145fm, obtained from the chirally extrapolated nucleon 
mass. For valence and sea quarks we use the clover action, with Cgw determined non-perturbatively. 
The charm quark mass was set by tuning the spin-averaged charmonium mass |(/«,^^ + 3my/>j/) to 
its experimental value. The pion on these lattices is rather heavy: niji = 1007(2) MeV. However, 
this reduces the mass gap between the rj and the rjc which might enhance mixing effects. Com- 
putations took place on the local QCDOC using the Chroma software library [^. We utilized 100 
effectively de-correlated configurations, where on each 200 independent stochastic estimates were 
created. 

3. All-to-all propagator techniques 

In order to calculate the quark-antiquark annihilation diagrams, we require propagators from 
any start to any end point. An exact inversion of the Dirac Operator M is not feasible in terms of 
the memory requirements and computation time. Hence we calculate unbiased stochastic estimates 
of these all-to-all propagators. We define random noise vectors \ri'), i = I, . . . ,N with components. 




v,w G {±1}. 



(3.1) 



We solve the linear problem for these N sources: 



s') =M-'\'n') 



(3.2) 
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where, 



l^L<aM%- = S.,y8a,b8,,p + ol^-^^ ■ (3.3) 



The propagator M ^ can now be estimated as follows: 



(3.4) 



We obtain the full propagator plus off-diagonal noise term, which vanish like We apply 
several improvement schemes to reduce this stochastic noise. 

3.1 Staggered spin dilution (SSD) 

Usually a propagator receives its dominant noise contributions from the stochastic sources 
that are closest to the sink, in terms of Euchdean distance or in spinor space. Such contributions 
can be reduced in a straight forward way by partitioning (or "diluting") spacetime, colour and/or 
spin into np disjoint subspaces |Jl^ |ll[]. The subsequent reconstruction of an all-to-all propagator 
from partition-to-all propagators comes at a significant computational overhead (proportional to 

-1/2 

np). Sometimes the resulting statistical error however is reduced by more than a factor ex np , 
justifying this method. 

Spin partitioning is usually achieved by employing four sets of noise vectors, each with only 
one spin component different from zero. Obviously, this is not the only possible pattern. We find 
that especially for heavy quarks other distributions of the spin components are superior to this stan- 
dard scheme. We attempt to exploit the fact that for large quark masses the Wilson-Dirac operator 
only weakly couples the upper and the lower spinor components and we devise a partitioning pat- 
tern where at different spacetime positions different spin components are stochastically seeded: an 
upper component site should only neighbour lower component sites and vice-versa. We found a 
scheme to be particularly effective, which we coin Staggered spin dilution (SSD). This is defined 
and compared to the standard dilution scheme in the left column of Fig. |l]. The numbers indicate the 
spinor component which is different from zero. The right column of Fig. |l| shows the corresponding 
nearest neighbour coupling strength, where red means "strong" and green means "weak". 

3.2 Hopping parameter acceleration (HPA) 



This technique [12-14] is based on the Hopping Parameter Expansion (HPE) of a Wilson- 



type Dirac operator M = 1 — kD, where D only couples nearest neighbours: 

oo 

M"^ = (I-kD)-^ =l + KD + ... + {KD)"-^ + {KD)"Y,{^Dy (3.5) 

(=0 

=> {kD)"M-^ =M-^ -{\ + kD + ... + kD"-^) : (3.6) 



the multipication of the propagator by (kD)" subtracts all the terms up to « — 1 hops, cancelling 
the corresponding fluctuations from noisy source positions near the sink. While this reduces the 
variance, the signals of Green functions for distances larger than n — I hops remain unaffected. 
Applications are the calculation of correlation functions with time separations t >na or of discon- 
nected loops where the first n — I terms are either calculated exactly or vanish identically. The 
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Figure 1: The left column shows the dilution schemes in two dimensions. The corresponding coupling is 
sketched in the right column (green is weak, red is strong coupling). 



effect of HPA on the pseudoscalar loop calculated using the Wilson action at the charm mass is 
plotted on the left hand side of Fig. (The first multiplication has no effect since Tr/s is zero 
anyway.) This technique is particularly suitable for heavy quarks where the rate of convergence of 



HPE improves but significant gains have also been reported for light quarks [|15[]. 

3.3 Recursive noise subtraction (RNS) 

In addition we use an algebraic improvement scheme (RNS) where the off-diagonal terms 



appearing in Eq. ( |3.4D are calculated and subtracted by hand. The idea can be illustrated as follows: 



= |.)(T]| +M-Hl-|77)(T]|) ^ |5)(77| + |5)(77|(1-|TJ)(77|). (3.7) 



The outer product |?7)(?7| should be truncated to diagonal blocks of a dimension that is small 
compared to a//V. We use the 12-dimensional colour (g) spin subspace and display a scatter plot 
between estimates on the right hand side of Fig. ^ The correlation angle is close to 7r/4: adding the 
two terms as suggested will result in a reduced stochastic error. In principle one can make a more 

general ansatz like |5')(77| [1 + a(l — |s)(t]|) + j8(l — |5')(t] |)^ H ] and optimize the parameters 

a,j5,... accompanying the different estimates of zero to minimize the noise in a given channel. 
Here we use a = 1 , j3 = 0. 

3.4 Truncated solver method (TSM) 

For the computation of the light quark loops appearing in the mixing matrix we apply the TSM 
using the BiCGStab solver. For details see [[T^. 

3.5 Reduction of the total computational effort 

We summarize the effect of the tested improvement schemes on the disconnected part of the 
zero-momentum projected tJc two-point function ( Tr(M"^75) Tr(M" Vs) ) at ? = in Table [l[ The 
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Figure 2: The effect of the HPA on the error of the pseudoscalar correlator for the Wilson action at the 
charm mass is shown on the left. On the right the RNS scatter plot is displayed. 



gain is in terms of real computer time, n is the power of kD applied to the sink vector. Our use of 
the clover action to calculate Tr(M^'75) restricts us to n < 3, since Tr(D^75) oc FF ^ 0. In principle 
one could calculate these terms by hand and add them back in again. With colour and staggered 
spin dilution and two applications of kD we obtain a net gain factor of almost 12. Note that these 
numbers do not include the effect of the TSM that we only use for the light quark propagators. 
Also we have not yet combined RNS with the partitioning methods. 
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no 


spin 


color 


color + spin 


SSD 


SSD + color 


RNS 
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1.43 
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1.97 


3.63 


1.87 


2 


2.89 


6.32 


5.06 


10.24 


7.16 


11.80 


5.44 



Table 1: Effective gain of the different improvement schemes, n is the power of kD used in HPA. 



4. The mixing 
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Figure 3: The mixing matrix. The blue lines indicate charm quark, the red ones light quark propagators. 

Now that we are equipped with techniques to address quark-antiquark-annihilation diagrams 
reliably, we calculate the elements of the mixing matrix. Both charmonium and light meson in- 
terpolators are included, three in each sector: (cc)o, (cc) lo, (cc)%q, {uii)o, {uu)=,, {uu)4q, where the 
subscript indicates the number of Wuppertal smearing steps with 5 = 0.3, employing spatial APE 
smeared («ap£ = 15, a = 0.3) parallel transporters. These smearing parameters were determined 
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by optimizing the effective masses within the two sectors, ignoring disconnected contributions, see 
[§]. In Fig. m we sketch the two by two mixing matrix, where the blue lines represent charm and the 
red lines light quark propagation. The prefactors are due to the two mass degenerate light flavours. 

The variational method is applied to the mixing matrix by solving a generalized eigenvalue 
problem, C(fo)"^^^C(?)C(fo)"'/^ =^"it,to) Xj/", see [§] for details. For sufficiently large times 
the eigenvalues and -vectors will approach their asymptotic values. The components of a given 
eigenvector can be interpreted as the coupling strengths of the corresponding operators to the state 
under consideration. Their magnitude will provide us with information about the magnitude of 
mixing in the system. We proceed as follows: we first determine the eigenvalues of the 3 by 3 
submatrices separately within each of the flavour sectors, where for the moment being we ignore 
the disconnected contribution in the charmonium sector (see the left hand side of Fig. The light 
T] and its first radial excitation (77') are the lowest two eigenvalues of the submatrix containing 
only light interpolators and rjc and rj'^ within the charmonium sector. We find a diagonalisation of 
the full 6 by 6 matrix to be numerically unstable and hence restrict ourselves to the basis of the 
states (cc)io, {cc)^o, {uu)^ and (mm)4o for the full-fledged mixing analysis. In the right hand side 
of Fig. ^ the effective masses of the lowest three eigenvalues obtained from this basis are shown 
together with the eigenvalues obtained above, ignoring the mixing effects: within statistical errors 
no effect is seen. The eigenvectors contain more detailed information about the mixing. We display 
the components of the ground state T] eigenvector on the left hand side of Fig. |5| and those of the 
Tjc eigenvector on the right hand side. Indeed, the Tj does not contain any statistically significicant 
admixture from the 77^ sector and vice versa. 

We conclude that there is no significant T]-T]c-mixing and neither are there any other significant 
flavour singlet effects on the r]c mass. However, the unrealistically heavy pion mass might have 
affected our conclusion. Runs on 24^ x 48 lattices with w 400 MeV will clarify this issue. 
Furthermore, mixing with other states like glueballs deserves future attention. 




t/a 



Figure 4: The left plot shows the effective masses of the eigenvalues obtained from diagonalizing the 
corresponding submatrices. To the right these eigenvalues (in black) are plotted together with the ones 
obtained from the basis (cc)io, (cc)8o, (mm)5, (mm)4o. 
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Figure 5: The eigenvector components of the rj are displayed on the left, the ones of the rjc on the right. 
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